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ABSTRACT 

The  purpose  of  this  project  is  to  improve  the  ability  to  locate  infrasound  events  by  improving  the  Horizontal  Wind 
Model  (HWM)  used  in  the  prediction  of  infrasound  observables,  in  particular,  travel  time  and  back  azimuth.  Our 
general  approach  uses  new  wind  data  to  develop  a  higher-resolution  and  more  accurate  HWM  and  tests  the 
performance  of  Ground  to  Space  (G2S)  specifications  based  on  the  new  HWM,  relative  to  those  based  on  HWM-93 
using  a  large  set  of  ground-truth  events.  We  assembled  a  preliminary  data  set  of  180  events  with  which  to  exercise 
our  statistical  framework  for  testing  the  effect  of  atmospheric  wind  models  on  the  accuracy  of  infrasound  event 
locations.  We  also  assembled  wind  data  representing  a  100-fold  increase  compared  to  those  used  for  HWM-93  and 
incorporated  these  wind  data  into  the  mechanism  used  to  generate  improved  HWMs.  We  applied  this  mechanism  to 
construct  a  preliminary  version  of  HWM-07,  and  demonstrated  that  this  preliminary  model  displays  some  important 
features  that  were  missing  from  HWM-93,  such  as  wind  velocities  of  up  to  70  m/s  at  an  altitude  of  1 10  km.  We 
created  an  efficient  mechanism  to  trace  rays  through  atmospheric  models  for  a  large  suite  of  infrasound  source- 
station  pairs.  Initially,  we  have  based  our  generation  of  travel-time  and  back-azimuth  predictions  on  the  HARP  A 
3-D  ray-tracing  programs.  We  refined  the  statistical  framework  for  testing  the  effect  of  atmospheric  wind  models  on 
predictions  and  applied  it  to  first  arrivals  for  a  preliminary  data  set  of  251  source-receiver  pairs  from  108  events  to 
evaluate  and  compare  predictions  made  against  the  existing,  baseline  atmospheric  specifications  for  our  study, 

G2S  and  HWM-93/MSISE-00,  on  which  it  is  based,  above  50  km.  As  expected,  we  find  that  G2S  more  often  and 
more  accurately  predicts  stratospheric  first  arrival  times  and  back  azimuths,  whereas  the  two  models  are  both 
equally  poor  at  predicting  thermospheric  arrival  times  and  back  azimuths.  It  is  for  thermospheric  arrivals  that  we 
expect  the  most  improvement  to  result  from  the  improvements  of  HWM-07. 
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OBJECTIVES 


The  primary  purpose  of  this  project  is  to  improve  our  general  ability  to  locate  infrasound  events  by  improving  the 
HWM  used  in  the  prediction  of  infrasound  observables.  We  do  this  using  recent  wind  measurements  and  direct 
infrasound  observations. 

A  second  purpose  of  this  project  is  to  evaluate,  statistically,  the  abilities  of  new  and  existing  models  to  locate 
infrasound  events  and  predict  the  principal  infrasound  observables,  arrival  time  and  back  azimuth.  We  do  this  by 
comparing  model  predictions  against  a  large  set  of  infrasound  observations. 

RESEARCH  ACCOMPLISHED 


Background 

The  source  location  problem  is  a  central  element  of  nuclear  explosion  monitoring.  In  regional  scenarios,  it  can 
provide  important  corroboration  of  explosive  events  but  only  if  signals  can  be  reliably  attributed.  Source  location 
and,  thus,  location-related  observables  are  central  to  source  attribution.  In  remote  regions  that  are  not  well  monitored 
by  other  land-based  technologies,  such  as  the  southern  oceans,  infrasound  may  provide  the  only  information  that 
constrains  location. 

The  reliability  of  infrasound  source  locations  depends  on  the  reliability  with  which  signal  features  can  be  predicted 
from  source  properties.  This  requires  detailed  knowledge  of  sound  propagation  in  the  atmosphere  and  the  physical 
properties  that  control  it,  in  particular,  the  static  sound  speed  and  wind  field.  The  state  of  the  art  in  representing  the 
atmospheric  properties  pertinent  to  infrasound  propagation  is  the  Naval  Research  Laboratory  (NRL)  G2S  model 
(Drob  and  Picone,  2000;  Drob  et  al.,  2003).  This  model  combines  climatological  data  captured  in  MSISE-00  and 
HWM-93  with  numerical  weather  prediction  (NWP)  data,  such  as  the  Navy  Operational  Global  Atmospheric 
Prediction  System  (NOGAPS)  and  the  National  Oceanic  and  Atmospheric  Administration  (NOAA)  National  Center 
for  Environmental  Prediction  Global  Forecast  System  (NCEP-GFS).  Wind  models  that  are  based  solely  on 
climatology  often  substantially  underestimate  tropospheric  and  stratospheric  wind  jets,  which  leads  to  significant 
under  prediction  of  the  fraction  of  energy  that  is  ducted  through  these  layers  (Drob  and  Picone,  2000;  Drob  et  al., 
2003).  The  utilization  of  the  near-real-time  NWP  data  in  G2S  significantly  improves  estimates  of  tropospheric  and 
stratospheric  wind  jets  that  profoundly  influence  sound  speed  at  lower  altitudes  (<50  km). 

While  the  G2S  models  are  the  most  comprehensive  atmospheric  representation  available,  significant  infrasound 
residuals  for  thermospheric  arrivals  are  still  obtained  using  them.  Le  Pichon  et  al.  (2005)  found  that  back-azimuth 
residuals  using  G2S  models  averaged  several  degrees  for  repeated  thermospheric  propagation  paths.  Others 
(e.g.,  O’Brien  and  Shields,  2004)  observe  that,  where  only  thermospheric  arrivals  are  predicted,  the  predicted  times 
systematically  over-predict  observed  arrival  times.  The  principal  source  of  these  discrepancies  appears  to  be  the 
wind  model  in  the  mesospheric/lower  thermospheric  (MET)  region,  where  weather  data  cannot  assist.  Recent  direct 
measurements  of  winds  at  higher  elevations  indicate  that  HWM-93  significantly  underestimates  winds  there. 

Stronger  MET  winds  imply  larger  azimuthal  deflections,  lower  turning  heights,  and  shorter  travel  times  for  predicted 
arrivals. 


HWM-07  Prototype 


The  purpose  of  this  task  at  this  stage  of  the  project  is  to  demonstrate  that  we  can  successfully  incorporate  all  the  new 
atmospheric  wind  data  into  the  mechanism  to  be  used  to  generate  HWMs.  A  new,  unified  HWM  empirical  model 
framework,  applicable  to  all  altitude  regions  and  allowing  for  greater  extensibility  (e.g.,  future  increases  in 
resolution)  has  been  developed  and  is  in  intermediate  testing  phases.  The  revamped  HWM  horizontal  formulation 
now  includes  a  complete  linear  and  orthogonal  spatiotemporal  spectral  basis  set: 


5=0  /=0  n-l 
S  M  N 

5=0  m=0  n-l 


(1) 
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In  this  equation  U(%X,5,0)  represents  the  zonal  wind  field  at  a  given  model  level  as  a  function  of  day  of  the  year, 
longitude,  local  solar  time,  and  latitude,  respectively.  The  vector  spherical  harmonics  (VSH)  are  expressed  in  the 
complex  spectral  amplitudes  coefficients,  Crj,Br,i.  A  similar  relationship  for  V(%X,5,0)  exists  given  the  basis  function 
parity/orthogonality  relationship  U:  and  V:  {Bi.,Bi,-Ci-,-Ci} .  The  Fourier  time  modulated  VSH  basis 

functions  are 
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[cos(mA)  cos(^t)  +  cos(mA)  sin(^r)]  (5) 


where  Pn,m{^)  is  the  normalized  associated  Legendre  polynomials  and  (s,l,m,n)  are  the  seasonal,  local  time, 
longitudinal,  and  latitudinal  spectral  wave  numbers,  respectively.  For  maximum  wave  numbers  of  S  =  2,  L  =  3, 

M  =  2,  and  N  =  8,  there  are  580  non-degenerate  unknowns  per  model  level. 

The  vertical  variations  in  the  HWM-07  model  are  represented  by  5  km  B-spline  weighting  kernels  (e.g.,  De  Boor, 
1968),  as  illustrated  in  Figure  1.  The  HWM-07  empirical  model  specifications  for  the  zonal  and  meridional  wind 
components  are  thus  given  by  the  equation  [U,V]  =  llai(z)Ui(T,A,  8, 0),  where  ai(z)  is  the  amplitude  of  the  respective 
weighting  kernel  in  the  altitude  domain  di. 

The  new  model  formulation  and  implementation  also  trivializes  future  increases  in  spectral  resolution  and  input  data 
set  density  for  the  foreseeable  future.  Linearization  allows  the  replacement  of  the  non-linear  Levenberg-Marquard 
parameter  estimation  procedure  with  a  straightforward  linear  least-squares  estimation  procedure  (e.g.,  Rogers, 

1996): 


X  =  X,  +{K^S;^K+S-:rrS-;{y-KxJ 


-K-l  i^T  0-1/ 


(6) 


s^iK^s;'K+s;'y' . 

In  this  equation  X  is  the  model  parameter  estimates;  is  an  a  priori  estimate  of  the  model  spectral  parameters 


(7) 


[i.e.,  a(z)^Cr,  Ci,  Br,  based  on  a  previous  iteration 
and  data  subset;  variable  K  is  the  forward  spectral 
matrix  [i.e.,  the  model  basis  functions  ^t,X,8,0)\, 

S  is  the  covariance  matrix  of  the  model  parameter 
estimate;  Sa  is  the  covariance  matrix  of  the  a  priori 
model  parameter  estimate;  y  is  the  observations;  and  Sg 
is  the  covariance  matrix  or  inverse  uncertainty  of  the 
observations. 

The  data  used  to  generate  the  new  models  are  briefly 
summarized  in  O’Brien  et  al.  (2006).  The  entire  set  of 
O(IO^)  data  is  far  too  large  to  be  fit  at  once.  Instead,  we 
update  the  coefficients  iteratively,  at  each  step 
introducing  a  new  random  sub-sample  of  about  70,000 
to  160,000  data  from  the  entire  set,  depending  on  the 
altitude  region.  This  limits  working  memory 
requirements  to  less  than  10  GB  yet  allows  all  data  to 
eventually  influence  the  model  coefficients.  With  each 


B-spline  Basis  Functions 


Figure  1.  Example  of  B-spline  weighting  kernels. 
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iteration,  the  coefficient  variances  decrease;  that  is,  the  coefficients  improve,  so  that  their  influence  in  the  next 
iteration  increases  over  that  of  the  next  data  sub-sample.  Iteration  continues  until  the  coefficients  converge.  The 
number  of  data  sampled  from  each  data  set/type  and  altitude  region  is  made  roughly  equal  in  each  iteration.  Data 
sets  that  span  sparse  regions,  or  are  the  only  data  set  in  a  given  domain,  are  weighted  higher.  Conversely,  data  sets 
that  overlap  with  many  other  data  sets  are  weighted  lower. 

Codes  for  the  new  model  formulation,  including  a  new  client  module  that  will  be  distributed  to  users  when  the 
model  coefficients  are  more  mature,  are  more  compact,  succinct,  understandable,  efficient,  and  extensible  than  the 
legacy  HWM-93  code.  The  new  client  module  also  supports  OpenMP  parallelization  of  key  subroutines  and  do 
loops,  allowing  for  exploitation  of  recent  advances  in  low-cost,  multi-core  processor  architectures  and  supporting 
FORTRAN  compilers.  This  will  be  important  further  into  this  project  when  large  ensembles  of  infrasound 
propagation  calculations  for  ground-truth  events  are  required  to  refine  and  tune  model  parameters  and  performance. 

Preliminary  Model  Comparisons — HWM-07  versus  HWM-93 

Although  HWM-07  can  be  regarded  as  only  preliminary  at  present,  it  already  shows  markedly  improved  agreement 
with  the  newer  data  over  HWM-93.  As  anticipated,  the  preliminary  HWM-07  has  much  higher  amplitude  tidal 
winds  in  the  85  to  135  km  altitude  regions,  in  part  because  the  representation  of  the  diurnal  migrating  tide  is  now  in 
line  with  more-recent  measurements  and  scientific  understanding.  It  also  includes  the  observed  seasonal  variations, 
phase  quadratures,  and  hemispheric  symmetries  and  asymmetries  of  the  respective  zonal  and  meridional  wind 
components.  Typical  improvements  are  illustrated  in  the  top  and  bottom  rows  of  Figure  2,  which  shows  WINDII 
wind  measurements  (Shepherd  et  al.,  1993)  averaged  in  overlapping  2-hour  bins  and  10  degree  latitude  bins  (black 
lines),  together  with  corresponding  point-for-point  averages  from  HWM-93  (blue  lines)  and  the  current  preliminary 
HWM-07  (red  lines).  The  top  row  illustrates  the  strong  zonal  winds  observed  on  either  side  of  the  meridional  wind 
maxima  that  are  better  captured  in  HWM-07  but  are  almost  absent  in  HWM-93.  At  higher  latitudes,  some  winds 
now  captured  exceed  70  m/s.  However,  there  are  fitting  issues  still  to  be  resolved  in  HWM-07.  The  bottom  row  of 
plots  illustrates  the  strong  meridional  winds  observed  at  subtropical  latitudes  (±22.5°).  These  are  clearly  captured 
better  by  HWM-07  than  by  HWM-93,  but  not  fully,  thus  motivating  further  model  adjustments. 

In  addition,  HWM-07  contains  fewer  features  believed  to  be  artifacts  than  does  HWM-93.  For  example,  the 
amplitude  of  semidiurnal  tidal  oscillations,  in  particular  the  zonal  component  near  the  equator  above  115  km,  is 
significantly  weaker  in  HWM-07  than  in  HWM-93,  as  illustrated  in  Figure  3.  HWM-07  much  better  matches  the 
WINDII  satellite  data  in  the  left  panel,  which  show  no  evidence  of  the  strong  diurnal  component. 


December  March  June 


local  solar  time  [h] 


Figure  2.  Examples  of  WINDI  wind  data  (Shepherd  et  al.,  1993 — black  lines)  averaged  in  one-hour  time-of- 
day  bins  showing  significantly  improved  predictions  of  HWM-07  (red  lines)  over  HWM-93  (blue 
lines).  Top:  HWM-07  captures  most  of  zonal  winds  averaged  between  30S  and  40S  at  an  altitude  of 
105  km,  while  HWM-93  captures  little.  Bottom:  HWM-07  captures  more  of  the  meridional  winds 
between  lOS  and  20S  than  HWM-93,  but  further  improvement  is  still  possible. 
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WINDII 
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Figure  3.  Comparison  of  December  zonal  winds  [m/s]  at  130  km  altitude  from  WINDI  wind  data  (Shepherd 
et  al.,  1993)  with  predictions  of  HWM-93  and  HWM-07.  Spurious  tidal  harmonic  artifacts  that  are 
present  in  HWM-93  are  eliminated  in  HWM-07. 


Most  importantly,  more  accurate  wind  specifications  of  greater  magnitude  in  the  95  to  120  km  altitude  region  will 
result  in  significantly  lower  thermospheric  ducting  heights  for  infrasound.  This  will  lead  to  shorter  predicted  travel 
times  and  less  molecular  attenuation.  Travel  time  calculations  utilizing  current  HWM-93  specification  generally 
over  predict  thermosphere  wave  guide  travel  times  (e.g.,  O’Brien  and  Shields,  2004).  For  example,  arrival  times 
predicted  by  ray  tracing  though  G2S  for  four  of  the  five  observed  arrivals  from  the  Watusi  chemical  test  explosion 
on  September  9,  2002,  are  over  predicted  by  200  to  700  seconds;  the  first  ray  predicted  to  land  along  the  source- 
receiver  azimuth  of  the  fifth  station  does  so  about  60  km  past  the  station.  All  of  these  rays  are  predicted  to  turn  at 
heights  of  105  to  140  km.  In  addition,  given  the  current  frequency-dependent  infrasonic  attenuation  coefficients 
(e.g.,  Sutherland  and  Bass,  2005),  current  PE  and  ray-tracing  calculations  of  many  ground-truth  events  where 
thermosphere  arrivals  are  observed  predict  too  much  attenuation.  For  the  case  of  significant  cross-track  lower 
thermosphere  wind  conditions,  the  new  wind  specifications  will  result  in  the  calculation  of  greater  azimuthal 
deviations,  again  in  line  with  infrasonic  observations  (e.g.,  Le  Pichon  et  al.,  2006). 


Statistical  Assessment  Method 


The  virtues  and  shortcomings  of  any  wind  model  must  be  evaluated  with  respect  to  the  accuracy  with  which  it 
predicts  infrasound  data.  In  this  respect,  the  most  reliable  observables  are  arrival  time  and  back  azimuth.  O’Brien 
et  al.  (2006)  developed  statistics  for  quantifying  the  overall  abilities  of  models  to  predict  these  for  an  arbitrary  set  of 
source-receiver  combinations.  In  particular,  statistics  comparing  predictions  from  two  models  p  and  q  have  the  form 


(8) 


where  is  a  sum  of  squares  of  the  Nk  residuals  from  model  k\  Fp.q  can  be  built  from  any  subset  of  observations, 

such  as  travel  times  or  back  azimuths.  Later,  we  compute  Fhwm-g2s  for  HWM93-MSISE00  and  G2S  using  only 
travel  times,  only  back  azimuths,  and  both.  If,  as  we  believe,  G2S  better  predicts  the  observations,  F  will  be  greater 
than  1 .  Our  assertion  is  that  Fpq  is  approximately  distributed  as  ^  ,  against  which  we  can  test  its  statistical 

significance.  This  depends  on  having  good  estimates  for  the  residual  uncertainties.  We  estimate  these  for  fractional 
travel  time  T and  back  azimuth  deviation  from  the  receiver-source  azimuth^  using  the  following  error  models: 


(9) 

(10) 
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Here,  cr^,  and  (Ja  are  uncertainties  in  observed  arrival  time  and  back  azimuth;  <7r,  (Je,  and  (Jt  are  the  uncertainties  in 
source  time,  location,  and  elevation.  These  are  one  of  the  subjects  of  the  next  section.  The  actual  uncertainties 
introduced  into  A  and  Thy  the  model  are  Gmt  and  Gma,  which  we  seek  to  quantify.  We  assume  that  uncertainty  in 
station  location  is  relatively  insignificant. 

Preliminary  Infrasound  Event  Data  Set 

The  meaningfulness  of  assessment  statistics  grows  with  the  number  of  observations  that  can  be  compared  to  model- 
based  predictions.  Thus,  a  significant  component  of  this  project  is  devoted  to  collection  and  analysis  of  ground  truth 
event  data.  In  principle,  any  event  is  usable  within  our  framework,  as  long  we  can  reasonably  estimate  cr^,  o;,  and  (J„ 
which  strongly  influence  Ga  and  Gt  and,  hence,  the  statistical  power  of  each  datum  and  its  weight  in  the  assessment 
statistics.  Without  accurate  uncertainties,  we  risk  polluting  statistics  with  erroneous  data  whose  value  is  erroneously 
overestimated.  Thus,  a  significant  component  of  the  ground-truth  analysis  is  devoted  to  estimating  the  source 
parameter  uncertainties.  Typically,  the  less  well  we  can  estimate  the  source  parameters  themselves,  the  less  well  we 
can  estimate  their  uncertainties,  so  we  conservatively  eliminate  poorly  constrained  events,  such  as  large  earthquakes 
and  rocket  launches,  having  spatially  and  temporally  extensive  sources.  Other  events  with  large  estimated 
uncertainties  may  be  eliminated  subsequently  if  they  appear  to  degrade  the  statistics  disproportionately. 

The  map  in  Figure  4  shows  source-receiver  paths  for  the  180  events  for  which  we  have  estimated  cr^,  o;,  and  Gj  to  be 
small  enough  to  admit  either  travel  time  or  back  azimuth  to  the  observation  set.  Of  these,  108  are  in  western  North 
America  and  3 1  at  the  Zheleznogorsk  mining  site  in  Russia.  The  remaining  events  are  distributed  fairly  randomly 
around  the  globe. 

In  the  last  year  we  have  improved  the  ground-truth  data  for  over  120  events.  Source  parameters  for  all  the  events  are 
stored  in  the  Monitoring  Research  Program  infrasound  database  found  on  the  web  at  http://www.rdss.info. 
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Figure  4.  Map  of  events  in  the  preliminary  infrasound  event  data  set. 
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Atmospheric  Propagation  Modeling 

We  have  obtained  model-based  predictions  by  tracing  rays  through  the  candidate  models  using  the  HARP  A  3-D  ray¬ 
tracing  programs  packaged  with  InfraMAP  (Norris  and  Gibson,  2006).  We  have  not  attempted  here  to  account  for 
terrain  except  by  elevating  the  sources  and  receivers  to  their  actual  heights  above  sea  level.  Future  computations  will 
account  for  terrain  using  the  NRL  RAMPE  ray  tracer  and  the  terrain-correcting  versions  of  HARP  A  recently 
released  in  beta  with  InfraMAP. 

We  executed  eigenray  searches  for  108  events  for  which  we  have  already  generated  G2S  models,  including  640 
arrivals  from  25 1  source-receiver  pairs.  Even  for  this  fraction  of  the  anticipated  final  data  set  and  only  two 
atmospheric  models,  G2S  and  HWM93-MSISE00,  the  ray  tracing  represents  a  large  computational  task  -  it  required 
about  160  hours  of  CPU  time  on  an  AMD  Opteron  285  microprocessor,  although  instabilities  with  finding  the 
intersections  between  rays  and  the  ground  may  account  for  a  significant  portion  of  that  time.  The  task  also  demands 
significant  bookkeeping,  since  the  ray-tracing  programs  must  be  executed  separately  for  each  source-receiver-model 
scenario,  each  producing  multiple  outputs.  To  facilitate  the  repeated  task  of  computing  eigenrays  for  an  arbitrarily 
large  set  of  source-station  pairs,  we  have  created  efficient  UNIX  scripts  for  running  the  HARPA  programs, 
extracting  the  prediction  data  from  its  outputs,  and  assimilating  its  results  into  an  Oracle  database,  where  the 
evaluation  statistics  are  computed. 

However,  the  process  of  mapping  predictions  to  arrivals  is  more  difficult  to  automate  or  even  objectify.  In  some  way 
this  must  be  done  by  correlating  the  values  of  measured  and  predicted  observables,  such  as  arrival  time,  back 
azimuth,  slowness,  amplitude,  and  duration.  Even  some  unobservable  predicted  ray  properties  such  as  turning  height 
or  number  of  ground  reflections  may  influence  the  mapping.  It  is  not  always  clear  how  much  weight  each  predicted 
and  observed  datum  should  get  in  the  process,  creating  a  significant  combinatorial  problem,  especially  if  there  are 
multiple  observed  and  predicted  arrivals  for  a  given  source-receiver  pair. 

Ultimately,  we  will  have  to  face  these  vagaries,  but  here  we  present  statistical  comparisons  based  on  a  more 
illustrative  mapping  that  avoids  them  altogether.  We  map  only  first  arrivals,  each  to  the  predicted  arrival  with  the 
best  fitting  arrival  time.  The  mapping  does  not  guarantee  or  even  purport  to  be  the  best  mapping  in  the  sense  of  all 
the  observables — some  mappings  are  certain  to  be  incorrect — but  it  is  unambiguous.  We  also  consider  separate 
statistics  for  arrivals  predicted  to  have  maximum  turning  heights  in  the  troposphere/stratosphere  (TS)  region 
(<65  km)  and  the  MET  region  (>  65  km)  because  we  expect  predictions  based  on  G2S  to  differ  significantly  from 
those  based  on  HWM93-MSISE00  only  where  the  rays  spend  most  of  their  time  at  TS  altitudes.  Figure  5  (a) 
illustrates  that,  in  fact,  all  the  predicted  turning  height  are  easily  separated  into  these  height  regions,  all  having 
maximum  turning  heights  below  55  km  (TS)  or  above  100  km  (MET). 

As  expected,  using  G2S  leads  to  predictions  of  significantly  more  TS  arrivals  (about  1/3  more)  and  more  arrivals 
overall  than  HWM93-MSISE00.  This  is  illustrated  in  the  Figure  5(b)  histogram  of  the  percentage  of  source-receiver 
pairs  for  which  the  models  predicted  any  arrival,  grouped  by  the  ducting  region  predicted  for  the  first  arrival.  The 
anticipated  superiority  of  G2S  is  confirmed  in  Table  1,  in  which  we  have  tabulated  Fg2s/hwm  and  the  corresponding 
significance  levels  for  all  combinations  of  travel-time  and  back-azimuth  data  in  each  of  the  TS  and  MET  ducting 
regions.  From  all  perspectives  G2S  has  outperformed  HWM93-MSISE00  (Fg2s/hwm  >1)?  mostly  with  high 
statistical  significance.  The  histograms  of  travel-time  residuals  in  Figure  6  illustrate  that  a  major  reason  for  the 
domination  is  the  replacement  of  a  number  of  highly  erroneous  MET  first-arrival  predictions  of  HWM93-MSISE00 
by  much  more  accurate  TS  predictions  of  G2S.  This  leads  to  dramatically  better  predictions  in  both  ducting  regions, 
especially  of  travel  times. 


Table  1.  Fg2s-hwm  statistics  for  first  arrival  prediction  residuals  and  associated 
significance  levels  in  parentheses.  Values  greater  than  1  indicate  that 
G2S  residuals  are  smaller. 


Data  Type 

Maximum  Turning  Height 

<  65  km 

>  65  km 

All 

Travel  times  only 

1.84  (99.7%) 

1.20  (82.9%) 

1.40  (98.9%) 

Back  azimuths  only 

1.13  (67.5%) 

2.11  (99.9%) 

1.66  (100%) 

Travel  time  and  back  azimuths 

1.34  (95.7%) 

1.84  (100%) 

1.70  (100%) 
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Figure  5.  (a)  Turning  heights  of  predicted  first  -arrival  predictions  from  G2S  and  HWM93-MSISE00  define 
well-separated  troposphere/stratosphere  (TS)  and  mesosphere/lower  thermosphere  (MET)  turning 
height  domains,  (b)  Percentages  of  source-receiver  pairs  predicted  by  each  model  grouped  by  the 
turning  heights  of  the  predicted  first  arrivals.  G2S  predicts  significantly  more  lower-atmospheric 
arrivals.  We  anticipate  a  similar  increase  in  the  number  of  predicted  MET  arrivals  from  HWM-07. 


We  anticipate  that  HWM07-MSISE00  and  G2S-07  will  provide  additional  improvement  to  MET  predictions  over 
their  existing  counterparts.  Figure  6(b)  shows  the  strong  tendency  of  either  HWM-93 -based  model  to  over  predict 
travel  times.  Some  of  this  is  no  doubt  the  result  of  failure  to  predict  much  faster  TS  arrivals.  However,  some  is  the 
simple  result  of  over  predicting  the  MET  travel  times.  With  generally  stronger  MET  winds,  HWM-07  will  have 
higher  overall  effective  sound  speeds  in  the  MET  region.  This  will  have  a  two-fold  effect  to  shorten  MET  travel 
times  by  increasing  the  propagation  speed  and  lowering  predicted  MET  turning  heights,  thereby  shortening 
propagation  paths.  A  less  likely  improvement  is  that  increased  mesospheric  winds  will  effectively  raise  the  top 
stratospheric  duct  in  some  cases,  thereby  increasing  the  number  of  faster  stratospheric  predictions. 
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fractional  travel-time  residuals 
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fractional  travel-time  residuals 
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Figure  6.  Fractional  travel-time  residuals  for  first  arrivals  predicted  to  turn  at  heights  within  (a)  the  TS 

region  below  65  km  and  (b)  the  MLT  region  above  65  km.  Travel  times  in  both  regions  are  better 
predicted  by  G2S.  MLT  residuals  are  strongly  biased  toward  late  prediction,  which  HWM-07  is 
intended  to  address. 


CONCLUSIONS  AND  RECOMMENDATIONS 


We  have  begun  construction  of  the  new  atmospheric  wind  specification  HWM-07.  The  new  model  is  structurally 
and  computationally  superior  to  its  predecessor  HWM-93.  Already,  we  have  achieved  significantly  improved  fit  to 
many  wind  data  sets  gathered  since  HWM-93  was  created  and  reduced  model  artifacts.  Continued  refinement  of  the 
model  parameters  is  still  required  to  fit  all  the  salient  features  that  are  manifest  in  the  wind  data. 

We  have  refined  our  data  set  of  ground-truth  events,  carefully  assessing  source  parameters  and  their  accuracy.  We 
have  increased  the  number  of  constraining  events  to  almost  200  and  will  continue  to  expand  this  set.  We  have 
developed  statistical  and  computational  machinery  required  to  generate  arrival  predictions  by  ray  tracing  for  a  large 
set  of  source-receiver  pairs  to  test  against  the  observations.  We  have  applied  the  machinery  to  first  arrivals  from 
HWM93-MSISE00  and  G2S,  clearly  and  quantitatively  demonstrating  the  anecdotally  accepted  prediction 
superiority  of  G2S  arrivals. 
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